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Figure 2: Example of a perturbatively predicted two-particle energy spectrum for isospin 
(solid lines) with a value ## = 24 for the renormalized coupling and a mass ratio m a /m % = 2 
(with rriir = 0.23). The dotted lines show the free energy spectrum (3.9). The <r-resonance 
energy m a is marked by the horizontal line at W = 3m v . 
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Figure 3: Lines of constant mass ratios m a /m v and m v L in the k-J plane derived with the 
help of the scaling laws of ref. [22]. The solid line represents the ratio m a /m^ = 3. Our choice 
of parameters (k, J) near that line is symbolized by the thick dots. To the left and to the 
right, the region where elastic two-particle scattering is expected is bounded by dashed lines 
corresponding to m a /m v = 2 and m a /m v = 4, respectively. The dotted line corresponds to 
m^L = 3 for L = 16. 
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Figure 5: The upper (lower) curve shows a two-parameter fit to the inverse a (ir) propagator 
in momentum space. The data are for the 32 3 • 40-lattice at (k = 0.315, J = 0.01). The 
intercept with the abscissa yields — rh 2 a (— m\). 



34 




block mean value ( . . . ) b 



® 




DIAGONALISATION 




DETERMINATION OF THE 
ENERGY EIGEN VALUES 




Figure 6: Different strategies for determining the block mean values of the two-particle energy 
spectrum. 
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Figure 7: The self-adjusting fit of the (mean) eigenvalues A„(i,£ ) of the correlation function 
matrix V(t,t ) for isospin at (k = 0.315, J = 0.010, L = 32): The eigenvalues on the different 
time slices are connected by dotted lines according to the criterion of the highest possible 
collinearity of the eigenvectors. Crossing of lines implies a contradiction to the sorting of the 
eigenvalues with respect to their size, which is expected to be equivalent to the former method. 
The crossing point determines the upper end t max of the range used in the exponential fit (solid 
line). 



36 




Figure 8: Comparison of the Monte Carlo results (crosses) with the perturbatively predicted 
two-particle energy spectrum in the isospin-0 channel for the simulation point (k = 0.315, 
J = 0.01). The data are in excellent agreement with the perturbative predictions (solid lines) 
based on the results of table 11. Dotted lines refer to the free energy spectrum, while the 
dashed lines show the perturbative spectrum based on the estimates of table 4. The location 
of the resonance energy m a is symbolized by the dotted horizontal line at W^3m^. - Errors 
are smaller than the symbols. 
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Figure 9: Comparison of the Monte Carlo results (crosses) with the perturbatively predicted 
two-particle energy spectrum in the isospin-2-channel for the simulation point (k = 0.315, J = 
0.01). The data are in excellent agreement with the perturbative predictions (solid lines) based 
on results of table 11. Calculations based on the estimates of table 4 yield indistinguishable 
results. Dotted lines refer to the free energy spectrum. - Errors are smaller than the symbols. 
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Figure 10: All data for the two-particle energy spectrum for isospin (l.h.s.), and for isospin 
2 (r.h.s.). The dotted lines in the isospin-0 case show the free spectrum. For both isospins and 
all simulation points (k = 0.320, J = 0.020, top), (k = 0.315, J = 0.010, middle) and (k = 0.310, 
J = 0.005, bottom), the data are in excellent agreement with the perturbative calculations 
(solid lines) provided we use the resonance parameters of table 11 in the case of isospin 0. - 
Errors are smaller than the symbols. 
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Figure 11: Normalized spectral amplitudes \c v -\ 2 in the isospin-0 channel for the operators Oj, 
j = <j, 0, 1, . . . in the case of plane waves at (k = 0.315, J = 0.01, L = 24). 
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Figure 12: Normalized spectral amplitudes \c v -\ 2 in the isospin-0 channel for the operators Oj, 
j = <j, 0, 1, . . . in the case of Liischer's wave functions at (k = 0.315, J = 0.01, L = 24). 
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Figure 13: — <t>y^j mod7r as function of k/m v with L fixed for (k = 0.315, J = 0.01) and 
L — 16, 20, 24, 32. The function — (f> taken at the momentum k corresponding to one of the 
energy values of table 8 determines the scattering phase shift at this energy. The error of the 
scattering phase shift is calculated from the error of k using the slope of cf). The dashed lines 
represent the behaviour of the scattering phase shift according to a fit with respect to the 
perturbative ansatz (3.6) (compare with fig. 1). The dotted vertical line indicates the inelastic 
threshold. 
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Figure 14: Lattice effects in the calculation of the scattering phase shifts for isospin and 
2: the upper (lower) points are determined by means of the continuum (lattice) dispersion 
relation. The scattering phases in the isospin-0 channel are well explained by a fit with the 
perturbative ansatz (3.6), and for isospin 2 the data are in good agreement with perturbative 
calculations based on these results (solid curves). Here we show the results at the data point 
(k — 0.320, J = 0.02) with the largest lattice effects. 
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Figure 15: Comparison of the fit with respect to the perturbative ansatz (3.6) (solid curve) 
and the perturbative prediction based on the estimates rh a and gn of table 4 (dashed curve) 
in the isospin-0 channel for (k = 0.315, J = 0.01). 
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Figure 16: Scattering phase shifts for isospin 1 at (k = 0.315, J = 0.01). The data are in good 
agreement with perturbative calculations based on the estimates rh a and gn of table 4 (solid 
curves). 
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Figure 17: Scattering phase shifts for isospin (l.h.s.) and for isospin 2 (r.h.s.). In the case 
of isospin the data are very well explained by a fit according to the perturbative ansatz (3.6), 
and the values of the scattering phase shift in the isospin-2 channel are in good agreement with 
the perturbative calculations based on these results (solid lines). 
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1 Introduction 

Most experiments in High Energy Physics are scattering experiments, the analysis of 
which leads to scattering phases, especially in the elastic channel. Therefore, it is a 
major challenge for theory to calculate these phases (or the corresponding scattering 
amplitudes). In the past, this has been done extensively by means of perturbative 
methods. However, if couplings become strong, non-perturbative methods are required. 
The most prominent example of a strongly coupled theory is QCD at low energies: 
Up to now it has been impossible to derive hadron-hadron scattering phases in an ab 
initio calculation from the fundamental theory (except for pion scattering lengths in the 
quenched approximation on the lattice [1-4]); in particular, the formation of resonances 
like, e.g., the p-meson, is still poorly understood. 

The most promising non-perturbative approach in quantum field theory consists in 
putting the theory on a lattice and computing the functional integral by Monte Carlo 
simulations. Therefore, it is highly welcome to have a method that allows to extract 
elastic scattering phases from quantities which can be calculated in a Monte Carlo "ex- 
periment". Such a method has been devised by Liischer [5,6]. It is based on a relation 
between the energy spectrum of two-particle states in a finite box with periodic boundary 
conditions and elastic scattering phase shifts defined in infinite volume. This relation is 
proved to hold for the scattering of two identical bosons below the inelastic threshold 
in any massive quantum field theory provided polarization effects are negligible. Since 
the spectrum of two-particle states in a finite volume can be computed numerically, we 
are thus in a position to determine phase shifts from Monte Carlo simulations. The fact 
that numerical simulations are necessarily performed on finite lattices is no disadvantage 
here, rather the finite volume is exploited to probe the system. 

Scattering lengths have already been investigated [1-4,7-9] by means of the asympto- 
tic large- volume behaviour of the lowest two-particle energy levels [10,11]. Furthermore, 
Liischer's relation for two-dimensional models has been used to extract scattering phases 
in the 0(3) nonlinear cr-model [12] and to study resonance scattering for two coupled 
Ising systems [13]. It has also been applied in three dimensions to meson-meson scatte- 
ring in QED [14]. For previous attempts to investigate resonances on finite lattices see 
refs. [15,16]. 

In this paper we demonstrate the applicability of the procedure in the four-dimensional 
0(4)-symmetric $ 4 -model. The model is studied in its broken phase, and the partic- 
les whose scattering is investigated are the three Goldstone bosons ("pions"). Since 
Liischer's method only works in the absence of massless particles, we introduce an exter- 
nal source term which breaks the 0(4) symmetry explicitly and hence makes the pions 
massive. However, the pion mass is kept small enough to allow for a resonance (the "er- 
particle" ) in the appropriate channel and we are able to calculate the mass and the width 
of this resonance. Most likely our model is not strongly interacting in the scaling region, 
so we have perturbation theory to compare with. Indeed, we find impressive agreement. 
The discrepancy noted in our preliminary results [17] turned out to be caused by an 
inconsistent normalization [18]. 

The paper is organized as follows: Sect. 2 summarizes Liischer's method. The per- 
turbative results for the scattering phases are given in sect. 3. In sect. 4 we describe the 
details of our Monte Carlo simulation. After a discussion of the single-particle spectrum 
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in sect. 5 we turn to the calculation of the two-particle energies (sect. 6). The results 
of these calculations are presented in sect. 7, followed by a brief discussion of the wave 
functions used in the construction of the two-particle operators (sect. 8). In sect. 9 we 
apply Liischer's method and compute scattering phases from the measured two-particle 
energies. In particular, we obtain the mass and the width of the er-resonance. Sect. 10 
contains our conclusions. 



2 Theoretical background 

In this section we summarize the formulae relating the energy spectrum of two-particle 
states in a finite periodic box and elastic scattering phases. We shall not give this 
relation in full generality, but restrict ourselves to the cases which are relevant for the 
applications we have in mind. For more details as well as for the proofs we refer to 
ref. [5]. 

Consider a system of two identical bosons of mass m, > with zero total momentum 
in a cubic box of size L 3 . The two-particle states in the elastic region are characterized 
by centre-of-mass energies W or relative momenta k related by 

W = 2\jm\ + k 2 , A: = I*: |, (2.1) 

where W lies in the range 

2m v < W < 4m w (2.2) 

or, equivalently, with k satisfying 

< K V^m, . (2.3) 

Rotational invariance being broken by the cubic box, the two-particle states are classi- 
fied according to the irreducible representations of the cubic group 0(3, Z). Below the 
inelastic threshold, the discrete energy spectrum corresponding to a particular represen- 
tation is determined by the scattering phase shifts Si with angular momenta I which are 
allowed by the symmetry of the states. The subspace of cubically invariant states (i.e. 
states transforming according to the trivial representation Af) receives contributions 
from / = 0, 4, 6, . . . Similarly, the Tf -sector (vector representation) contains 1 = 1,3,... 
We shall only be interested in the cases I = and 1 = 1. Hence we have to study two- 
particle states transforming according to the representations A* and Tf of the cubic 
group. Fortunately the contributions of the higher angular momenta in each sector will 
turn out to be negligible. 

Assuming dominance of the lowest angular momentum in each symmetry sector, an 
energy value W v [y = 0, 1, 2, . . .) belongs to the two-particle spectrum in a finite volume, 
if the corresponding momentum k v = \J (W„/2) 2 — m\ is a solution of 



Si{K) = -<j ) [^- \ modvr (2.4) 
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for I = in the case of the ^-sector and for I = 1 in the case of the T x -sector (see 
ref. [6]). The function (f)(q) is defined as the continuous function satisfying 



for Re s sufficiently large and can be continued analytically to s = 1. The above 
statement holds only for (j^f) < 9, which is fulfilled in our simulations. Note that 
eq.(2.4) has been derived in the continuum, hence it may be applied to lattice results 
only if scaling violations can be ignored. Furthermore, the derivation neglects finite 
volume polarization effects (exchange of particles "around the world"). Whether this is 
justified, can be checked by a careful study of the single-particle states. 

In our Monte Carlo investigation we calculate the two-particle energy spectrum W v 
in a sector of definite cubic symmetry (Af or T-f) for a given size L. The scattering 
phase at the corresponding momentum can then be read off from eq.(2.4). In order to 
scan the momentum dependence of the phase shift we vary the spatial extent L of the 
lattice. 

Conversely, if 6i(k) is considered to be known (e.g. from perturbation theory), eq. (2.4) 
allows us to compute the momenta k v and hence the two-particle energy spectrum for 
a given box size L. Fig. 1 illustrates the graphical solution of eq. (2.4) using the lowest 
order perturbative result for S (k) (in the isospin-0 channel, see following section) as 
input: The momenta k v corresponding to the energy values W v are determined by the 
intersections of the solid curve ((f>(q)) with the dashed lines (6o(k) — v% for v = 0, 1, 2, . . .). 
Note that in this example the lowest energy state (v = 0) is calculated from the analytical 
continuation of the scattering phase to imaginary values of k. 

3 Perturbative results for the scattering phases 

Most likely, the $ 4 model in four dimensions is trivial, and hence the renormalized 
coupling is bounded from above (for a recent review see ref. [19]). Indeed the bound 
turns out to be so low that the theory can be considered as weakly interacting in the 
whole scaling region (at least in the lattice regularization used in this work) [20] . There- 
fore perturbation theory should be applicable and lead to predictions for the scattering 
phases. 

Considering the 0(7V)-symmetric $ 4 -theory we write the action in the euclidean 
continuum as 




#0) = 0. 



(2.5) 



The zeta function Z 00 (s; q 2 ) is given by (also for q 2 < 0) 




(2.6) 



S{$ 



;j,™ 2 } = J 



N 



(3.1) 



with a = 1,...,N. Since we work in the broken phase, the mass parameter m 2 is 
negative. The external field, which breaks the symmetry down to O(N-l), is taken to 
point in the TV-direction and has the magnitude j. 
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The pions, whose scattering is to be investigated, are characterized by their momenta 
k and an "isospin" index a = 1, . . . , N—l. The single-pion states transform according to 
the vector representation (isospin 1 for N = 4) of the residual O(N-l) symmetry. The 
product of two of these representations decomposes into three irreducible representations 
(with isospin I = 0,1,2, respectively, for JV = 4). The corresponding projectors Q 1 are 
given by 

Qa'b',ab = N _ 1 Sa'b'Sab , 

Ql'b',ab = 2 i^a'aSb'b ~ 8 a 'b$b'a ) (3.2) 
Ql'b',ab = 2 (fia'a^Vb + ^a'bh'a) ~ jy _ ^ a'b'^ab ■ 

With their help, the amplitude T for elastic two-pion scattering can be expressed in 
terms of the scattering amplitudes T 1 for fixed isospin 7: 

2 

t = y.Qi tI - ( 3 - 3 ) 

7=0 



In the centre-of-mass system with total energy W = 2ym£ + k 2 the amplitude T 1 is a 
function of the scattering angle 6 and the absolute value of the momentum k = \k\. The 
partial wave decomposition reads: 

1 Pi \Al 00 

T 1 = —r—Y,( 21 + 1) P <(costf) t\{k) , (3.4) 
ft 1=0 

where P; denotes the Legendre polynomials. Note that due to Bose symmetry the partial 
wave amplitudes t\ vanish, if 7 + / is odd. In the elastic region 2m ff < W < 4m ff we can 
express t\{k) in terms of the real scattering phase Sf(k): 

tj = ^(e^-l) (3.5) 

(provided the S-matrix is unitary). 

We are interested only in the case where the er-particle is unstable (m CT > 2m 7r ) and 
appears as a resonance in the corresponding channel (7 = 0, Z = 0). The scattering 
amplitude becomes singular and the perturbative calculation requires a separation of 
the singular contribution 8® s to 5$ from the regular piece ^ r [21]. 

In lowest nontrivial order one obtains the following results [21]: 

K = C + S l (3-6) 

co N - 1 ml - ml k 

with tantf 0i . = * fl — ^^-^j-, 

co c2 N — l k 

c 2 _ 9b_ rrij -mj / 4fc 2 + m 2 a \ _ 9b_ k_ ( ~ 7 \ 
d ° 967T kW V m 2 ) 24nW i 

ci 9R ml-ml ( ml, , 4k 2 + m\ , g R m\ - mj 

dl ~ 96tt kW y L+ 2k*r n y ml ) 48tt kW ' 
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Here m, and m a denote the physical masses of the particles (the real part of the propaga- 
tor pole in the case of the a) and is the renormalized coupling constant. Comparison 
of <5q, Si with Sf, Sf for / > 4 shows that the higher angular momenta are completely 
negligible in perturbation theory. The same statement holds in the 7=1 channel: 6} 
with I > 3 is much smaller than S\. 

The perturbatively calculated scattering phase 6° may be continued analytically to 
purely imaginary values of k. In the 7 = channel this leads to scattering states which 
are obtained as solutions of eq.(2.4) for negative q 2 = (kL/2Tr) 2 . 

Fig. 2 shows an example of a perturbatively predicted two-particle energy spectrum. 
The dotted lines indicate the free energy spectrum 



W 

— = 2, 



k 2 - 2tt 

5*1 + —, k = —n, n£Z 3 . (3.9) 
\ ml L 

4 The Monte Carlo simulation 

On the lattice we parameterize the action as follows: 

S{§; J, k} = ^2 { ~ K ( $ a (x)$ a (x + A) + $ a (x)$ a (x - ft) ) 

x I ii=l 

+ A ( <S> a (x)$ a (x) - 1 ) 2 + <S> a (x)$ a (x) - J$ N (x)^ , (4.1) 

Here we assume a hypercubic lattice whose spacing is put equal to 1. As usual, fx denotes 
the unit vector pointing in //-direction. The relation between the lattice parameterization 
(4.1) and the continuum parameterization (3.1) of the action is given by the formulae: 

r(x) = V2^$ a (x), j = -Lj, 

, 2A - 1 a 
-m 2 = h 8 , 



(4.2) 

2A - 1 6A 

~2 ' 



In the actual simulation we used N = 4 and A = oo so that the scalar field is represented 
as a real four-component vector $(a;) of unit length: $ a (x)$ a (x) = 1. 

We shall work on lattices of size L 3 ■ T using a cluster algorithm. The simulation 
parameters k, J, L, T have to be chosen judiciously in order to extract a maximum of 
information from a given amount of computer time. As guideline we use the results on 
the scaling behaviour of our model obtained in ref. [22]. 



4.1 Choice of the simulation parameters 

On the basis of the scaling laws (see ref. [22]) we plot in fig. 3 lines of constant mass 
ratios m a /m v in the (k, J)-plane. In order to be sure that the er-resonance lies in the 
elastic region 2m v < W < 4m„., we want to choose the parameters k and J such that 
m CT /m ff ~3 (solid line in fig. 3). Furthermore, we have to suppress vacuum polarisation 
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effects by making m v L sufficiently large. The dotted line in fig. 3 serves to illustrate this 
constraint: Above this line we have m v L > 3 for L = 16. Finally we ended up with the 
choices («, J) = (0.310, 0.005), (0.315, 0.010) and (0.320, 0.020) shown as black dots in 
fig. 3. 

The spatial extent L of the lattice should now be taken such that the corresponding 
momenta k (see eq. (2.4)) cover the interesting region around the resonance where the 
phase shift is about ir/2. For an estimate we put the resonance energy W = m a leading 
to a relative momentum 



I ml 



K = ^-m\. (4.3) 
Then the scattering phase shift b\ takes the value 7r/2, if 

L = ^ (4.4) 

with Q(o) such that (j)(q^) = n/2 (mod7r). According to eq.(2.5) these values are just 
the zeros of the zeta-function Z 00 (l; q 2 ). Table 1 contains a list of the first of these zeros. 



index 


n 2 
1(0) 


-1 


-0.095900719 





0.472894247 


1 


1.441591313 


2 


2.627007612 


3 


3.536619947 


4 


4.251705973 


5 


5.537700774 


6 


7.196263202 


7 


8.287953655 


8 


9.534531427 



Table 1: Zeros of the zeta-function Z 00 (l;q 2 ). 



(k, J) m a 7h a /m v L 



1 


(0.310, 0.005) 


0.5383(11) 


0.1844(9) 


2.88 


22, 


38 


2 


(0.315, 0.010) 


0.7202(13) 


0.2309(8) 


3.13 


16, 


27, 37 


3 


(0.320, 0.020) 


0.9059(16) 


0.2977(10) 


3.07 


13, 


22, 30, 35, 38 



Table 2: Numerically realizable lattice extensions L for which the two-particle energy spectrum 
has an energy level in the resonance region. 
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(k, J) LATTICE SIZES L 3 ■ T 



1 


(0.310, 0.005) 


16 3 


■32, 


20 3 


•32, 


24 3 


• 32, 32 3 • 40 


2 


(0.315, 0.010) 


16 3 


• 32, 


20 3 


• 32, 


24 3 • 32, 


24 3 • 60, 32 3 • 40 


3 


(0.320, 0.020) 12 3 • 32, 


16 3 


•32, 


20 3 


•32, 


24 3 


• 32 



Table 3: Actually realized simulation parameters (k, J) and lattice sizes L 3 ■ T. 



From test runs or from the scaling laws [22] we get estimates for m T , m a and hence 
k a . (Actually, we took for m a the mass m a extracted from a fit of the er-propagator in 
momentum space.) With the help of table 1, eq.(4.4) then leads to a list of L values for 
which a two-particle energy level close to the resonance energy is to be expected (see 
table 2). In particular, this list gives us the minimal useful L- values. Also taking into 
account that the computing resources limit L from above, we adjusted the simulation 
parameters such that at least four different lattice sizes could be used. Finally choosing 
T > L (to allow for a reliable determination of energy levels) we arrive at the simulation 
parameters shown in table 3. 

4.2 Technical details 

Our configurations are generated by means of the cluster algorithm [24,23] generalized for 
actions with finite external source: Once the multi-cluster structure {cj} is determined, 
in the standard multi-cluster method (for J = 0) the spins of every second cluster are 
reflected with respect to the hyperplane perpendicular to a randomly chosen global 
0(4)-direction r. Generalized to a finite external source J the clusters Cj are flipped 
with probability 

w($, Cj ) = ( l + exp(2Jr 4 £ (r ■ $(x)) ) ) \ (4.5) 

depending on the field configuration within the cluster Cj. 

The simulations are performed on lattices of size L 3 • T. On the HLRZ-CRAY 
YMP8/832 we can only run lattices of sizes 12 3 • 32, 16 3 • 32 and 20 3 • 32. The lar- 
ger lattices 24 3 • 32, 24 3 • 60 and 32 3 • 40 are simulated on the Landesvektorrechner 
Fujitsu/Siemens-Nixdorf SNI S600/20 in Aachen. The speed of our simulation program 
is about 180-200 MFlops on the CRAY and 500-650 MFlops on the Fujitsu. We got 
a mean computing time of about 6.0 • 10~ 6 seconds per lattice site and iteration on the 
CRAY and 1.8 • 10" 6 seconds on the SNI S600/20. 

As a rule we performed 150 000 to 330 000 iterations per simulation point (k, J) 
and lattice size L 3 ■ T, measuring after every fourth iteration. Averages over blocks of 
1024 measurements are stored on disk. Of course, a later increase of the block length 
is possible. In order to generate the random number sequences we used the Kirkpat- 
rick-Stoll-Greenwood (shift register) random number generator which vectorizes on both 
machines [25]. 
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Simulation points 





1 


2 


3 


K 


0.310 


0.315 


0.320 


J 


0.005 


0.010 


0.020 


m v 


0.1870(6) 


0.2300(6) 


0.2947(6) 


rho 


0.5383(11) 


0.7202(13) 


0.9059(16) 


S 


0.18455(3) 


0.23574(2) 


0.28096(2) 


z 


0.9745(5) 


0.9723(5) 


0.9715(5) 




0.9323(10) 


0.9130(12) 


0.8966(13) 


9r 


21.96(13) 


24.43(12) 


27.03(14) 


Cp 


-0.0614(4) 


-0.0634(4) 


-0.0602(2) 



Table 4: Infinite volume results used in order to compute the perturbative predictions of the 
two-particle energy spectrum. 

5 Pion mass and single-particle spectrum 

Defining the particle mass by the single-particle energy with spatial momentum p = 
in infinite volume we expect the relativistic energy momentum relation 

E(p) = yfml+p (5.1) 

also to be valid on a finite lattice provided m v L ^> 1 and \p\ <C 1. The first inequality 
insures the suppression of finite volume effects [26] and the second one is related to the 
lattice effects, which for our choice of lattice action decrease with the square of the lattice 
spacing. 

In order to control both polarization effects and scaling violations, we first examine 
the single-particle energies extracted from the exponential decay of the propagator 

( E E &(-n, t) ¥{n , 0) \ ~ e- E ^ (5.2) 

with p 2 = (2ir/L) 2 n 2 = (2ir/L) 2 j. Here $> a (n,t) is the spatial Fourier transform of the 
pion field Q a (x,t). For a general lattice function F(x) it is defined by: 

= 73 E e i2 * g - a ' L F(x). (5.3) 
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3 


n 


\P \/ m * 




E(p)/E lat (p) 


1 


(1,0,0) 


0.8539 


0.9932(12) 


0.9976(12) 


2 


(1,1,0) 


1.2076 


0.9924(8) 


0.9987(8) 


3 


(1,1,1) 


1.4790 


0.9893(8) 


0.9973(8) 


4 


(2,0,0) 


1.7078 


0.9885(6) 


1.0016(6) 


5 


(2,1,0) 


1.9093 


0.9859(4) 


1.0000(4) 



Table 5: Single-pion energy levels and their deviation from the relativistic continuum and 
lattice dispersion relation as function of \p \ for the simulation parameters (k = 0.315, J = 0.01, 
L = 32). In the second column we give the lattice momenta n with n 2 = j for which the energy 
is calculated. These momenta are used to compute p in eq. (5.5). 

The numerical data of the single-particle spectrum are well approximated by a fit of 
the form 

E(p) « y/ml + p* (l + c p p 2 ) (5.4) 
(see fig. 4). Alternatively the data can be compared with the lattice dispersion relation 

E(p) « E lat (p) = 2 asinh ( X -^ m \ + f ) , p, = 2 sin | . (5.5) 

A fit of the form (5.4) or (5.5) for fixed lattice extent L yields m ff (L) (and c p (L)). The 
L-dependence of c p (L) turns out to be so weak that we simply average over the results 
for different values of L to obtain our final estimate for c p in the infinite volume limit (see 
table 4). The pion mass, on the other hand, is extrapolated to infinite volume assuming 
a decay of the polarisation effects proportional to exp(— m 7r L)/(m 7r L) 2 . Both fits lead 
to consistent results for the pion mass. Table 4 contains the values obtained by means 
of (5.5). Determining by a fit to the single-particle energies we gain the advantage 
that all measured values of E(p) contribute to the calculation. 

Fig. 4 compares the measured single-pion energies E(p) at simulation point 2 (k = 
0.315, J = 0.01) with the lattice dispersion relation (5.5) and the continuum dispersion 
relation (some examples are listed in table 5). In both cases the infinite volume mass is 
used. Hence the plot demonstrates not only the quality of the fits (5.4) and (5.5) but 
also the smallness of the finite-size effects. 

Alternatively, the pion mass m v can be measured by a fit to the inverse propagator 
in four-dimensional momentum space: 

[G aa (p, -p)]- 1 = Z- 1 { ml + f } a = 1, 2, 3 . (5.6) 

The excellent quality of that fit is shown in fig. 5. After extrapolating to infinite volume 
we obtain within 1% deviation the same m^-values as with the determination from the 
single-particle spectrum. We also get the pion wave function renormalization constant Z . 
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For the perturbative calculation we need not only m v but also the er-mass m a and 
the renormalized coupling g^. The latter is defined by 



ml — m 2 



g R = 3Z °— \ (5.7) 

where £ is the field expectation value in infinite volume. We extract £ from the volume 
dependence of the expectation value of the magnetisation parallel to the external source 
as outlined in ref. [22]. The er-mass will later be determined from the analysis of the 
scattering phase shift 8®. However, a first estimate rh a can be obtained by a fit to 
the er-propagator in momentum space (see fig. 5). A volume dependence of m a is not 
observed. The results for m a as well as the ensuing estimates §r for the renormalized 
coupling are given in table 4, which also contains the wave function renormalization 
constant Z a extracted from the er-propagator. 



6 How to calculate the two-particle energy spectrum 

Following ref. [12] we extract the two-particle energy spectrum from a matrix of corre- 
lation functions of suitable operators, which couple to the two-particle states (with the 
appropriate quantum numbers). The states and operators are classified according to 
0(3)-isospin taking the values 0, 1 and 2 (see eq. (3.2)) and irreducible representations 
of the cubic group 0(3, Z) (Af or T-f). By definition the two-particle states are energy 
eigenstates and in general not identical to the states generated by the application of the 
operators to the vacuum. The spectral decomposition of the correlation function matrix 
will finally enable us to calculate the two-particle energy spectrum, if the number of 
operators exceeds the number of energy eigenstates below the inelastic threshold in each 
symmetry sector. 

6.1 Two-particle operators 

Operators for two pions with total momentum zero are defined by the following double 
sum over the spatial part of the pion fields [12] 

OfW-jjS E fi(x-y)^(x,t)^(y,t). (6.1) 

with suitable wave functions fi(x — y) (i = 1, 2, . . .). In order to save computer time we 
perform the actual calculations on the momentum lattice writing 

Of(t) = E kn) i> a (-n,i) &(n,t) (6.2) 

with the spatial Fourier transform defined in eq. (5.3). Since the scalar fields & a (x, t) are 
real, we have $ a (n, t) = $ a (— n, t), where a bar over a quanitity means complex conjuga- 
tion. This implies: Of(t) = 0\ a (t) for real wave functions fi(n) in momentum space 2 . 

2 0r fi(x) = fi(—x) in coordinate space. 
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The Fourier transformation <J? a (n , t) of the field is calculated by a three-dimensional Fast 
Fourier Transformation (FFT) which has excellent vectorization properties [27]. 

For isospin and 2 we compare two different kinds of cubically invariant wave func- 
tions. First we take plane waves (as refs. [12,13] in their two-dimensional models): 

fi(n) 

The index i = 0,1,2,... counts the integers j for which there is a three-dimensional 
integer vector n with n 2 =j. 3 

Our second choice is motivated by the periodic singular solutions of the Helmholtz- 
equation (see ref. [5]): 

(tf2_ 2\-l 

m= £ (fga ^)-i . * = 0,1,2,..., (6.4) 

n& 3 L 

where ^ = (L/2%)^(Wi/2) 2 — m 2 using as input the plane wave energy spectrum Wj of 
some test runs. Another way to calculate estimates ^ uses the results of table 4 together 
with perturbation theory: By computing the scattering phase shifts (3.6) and (3.7) and 
solving eq. (2.4) for ki one could estimate some values = kiL/2%. The wave functions 
(6.4) could have an enhanced projection especially in the isospin-0 channel at energies 
close to the resonance mass. 

Later we will show that the calculated two-particle energy spectrum does not depend 
significantly on the choice of the wave function. Hence both sets of operators allow a 
reliable determination of the energy values below the inelastic threshold. 

In the case of isospin 1 the analog of the cubically invariant wave functions (6.3) and 
(6.4) of the Ai sector are the vector-like wave functions ~ n Sj,n 2 or ~ n/(n 2 — q 2 ) of 
the T-f sector. For our calculations we have adopted the first possibility. 



S j,n 2 

E Sj,r 
mGZ| 



(6.3) 



6.2 Correlation function matrix 

For each isospin we examine the (connected) correlation function matrix of the two- 
particle operators 0^ h (t): 

c&t) ■■= Qi.», ab (ofW)of(o)) c 

= ( (WW)Of(0) ) - ( OfW) ) ( Of (0) ) ) . (6.5) 

Since we expect to find the er-resonance in the isospin-0 channel, we first restrict 
ourselves to that case. The corresponding two-pion operators 0{(t) are given by the 
trace of the operators (6.2) in the space of the internal 0(3)-symmetry (see eq. (3.2)): 

Oi(t) = ^=i:or(t). (6.6) 



3 j ^ i only for j > 7. 
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We also have to take into account the a field at zero momentum, since it has the same 
quantum numbers as the operators O; and is expected to create a state with energy 
below the inelastic threshold: 

o CT (t) = l> 4 (M = 7^ E * 4 (^M)- (6-7) 

With the help of eq. (3.2) we obtain from the definition of the correlation function 
matrix (6.5) the following expressions for the different isospin channels: 

Cijit) = (O i (t)O j (0)) c i,j = <t,0,1,2,..., 

Cl(t) = (lmOf(t)lmOf(0)) c i,j = 1,2,..., (6.8) 

Cg(f) = (ReOf(t)ReOf(0)) c -C^(t) i, j = 0,1,2, ... , 

where sums over repeated indices a, b are implied. Disconnected contributions are of 
course only expected for isospin 0. 

Instead of working with the connected correlation functions according to their de- 
finition we eliminate the disconnected contributions by measuring correlation function 
matrices of the following type 4 [12]: 

C ij {t) = ({O i {t)-O i {t + l)) 0,(0)} . (6.9) 

From the transfer matrix formalism one gets the spectral decomposition of these corre- 
lation function matrices: 

oo , 

Cijit) = J2Wvj e~ tW " with v] = VI - e~ w » ( v \Oj{0)\ Q ) . (6.10) 

The energy eigenvalues W v are assumed to be non-degenerate (in the symmetry sector 
under consideration) and are ordered such that W v < W„ +1 . The amplitudes v*f in the 
spectral decomposition (6.10) are proportional to the projections of the states generated 
by the two-particle operator out of the vacuum | Q ) onto the energy eigenstates \ v). 

Below the inelastic threshold 4m v there is only a finite number % of energy eigenstates 
W v . Furthermore, we can consider the correlation function matrix (6.10) only for a finite 
number f of indices This restriction is not only due to the finiteness of computing 
resources, but is also enforced by the fact that the number of linearly independent two- 
particle operators is finite on a finite lattice. In order to guarantee that the f-component 
vectors v u , v = 0,1,2, ... ,i — 1, are linearly independent f has to be larger than the 
number % of (expected) states below the inelastic threshold. 

For large time distances t the eigenvalues of the correlation function matrix C(t) 
are proportional to e~ tWv . Hence one could determine the energy eigenvalues from the 
eigenvalues of C{t) for t — > oo. However, since the statistical errors of the matrix C(t) 
are increasing with t, we use the following method to extract the energy spectrum. It 



4 Isospin indices are not shown explicitly in the following. 
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allows a reliable determination of the energy levels W v already for smaller values of t [12]. 
Given a (small) reference time t we consider the generalized eigenvalue problem 

C(t)w» = \ v (t,t )C(t )w v . (6.11) 

The eigenvalues X u (t,to) are given by 

A„(Mo) = e~^ w " (6.12) 

up to (negligible) corrections of order e~^ t ~ t ° >Wf . In order to analyse the generalized 
eigenvalue problem (6.11) we study the matrix 

V(t,t ) = C"5(i ) C(t) <H(i ) (6.13) 

having the same eigenvalues X v (t,t ) but eigenvectors u 1 " = (t )w". To calculate the 
spectral amplitudes we solve the equations 

Y,v? w? = 8^ e w " to l 2 (6.14) 

i 

by matrix inversion. Like (6.12) these relations suffer from exponentially suppressed 
corrections. 

The projection properties of the chosen operators onto the energy eigenstates are 
characterized by the normalized amplitudes 

< = / Vj > ^ = / 1 w < • ( 6 - 15 ) 



e 



They do not depend on the normalization of the operators Oj. 

For the calculation of the correlation functions we have also applied the two-cluster 
method [12]. But eventually we did not adopt this method because of two reasons: It 
takes too much computer time for other wave functions than plane waves (6.3) and for 
the most interesting isospin-0 case we do not expect (and actually did not find) any 
reduced variance [12]. However, we tested this method and obtained agreement with the 
standard measurements within statistical errors. 



6.3 Various methods of analysis 

In order to control systematic errors and to improve the quality of the analysis we ex- 
amined different methods to determine the two-particle energy spectrum (for a summary 
see fig. 6). They all lead essentially to the same results as is demonstrated in table 6 for 
a typical case. 

As described above, we diagonalize the matrix T>(t, t ) = (t )C(t)C~^ (t ) for some 
small value of the reference time t . This reference time has to be kept small in order to 
guarantee the numerical stability of the determination of the eigenvalues. 

To determine the energy spectrum W v we fit the eigenvalues with the exponential 
ansatz 5 

A„(Mo) = e- ( '-' o)W ". (6.16) 
5 Using a fit with the function cosh ( (t — to — T/2) W v ) does not affect the results. 
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The data are combined into blocks, we take the block mean values and calculate the 
total mean values and errors using the Jackknife method. The data are written out by 
the simulation program with a standard block length l^k of 1024 measurements. In order 
to get a correct estimate of the errors we increase the block length l U k until the errors 
stabilize. Already with l U k = 1024 the analysis seems to be stable, and we decided to put 
Zi/fc = 2048 in our final results. 

There are (at least) three different routes to take the mean of the data: We can 
combine into blocks the correlation function matrices T>(t,t ), their eigenvalues \ u (t,t ) 
or the energy levels W v from the exponential fit. These different strategies summarized 
in fig. 6 are denoted by (I), (2), and (3). A change of strategy only leads to negligible 
fluctuations of the energy levels (see table 6). Since the error analysis in strategy (3) 
seems to be most reliable, we use this method for calculating our final results. 

Another free parameter of the analysis is the dimension f of the correlation function 
matrix T>(t,t ), which has to be diagonalized. As long as f is larger than the number 
i of energy levels expected below the inelastic threshold we do not see any significant 
dependence of the energy spectrum on f. 

In order to check the dependence on the temporal lattice extent T, we have performed 
runs on a 24 3 • 60 lattice. The results do not differ significantly from those obtained on 
a 24 3 • 32 lattice. 

From the energy eigenvalues W v we have to calculate the corresponding momenta k v . 
Neglecting lattice effects one could use the continuum dispersion relation: (W v /2) 2 = 
m^ + k 2 ,. Alternatively, as for the single-particle spectrum, we can try to suppress the 
dependence on the lattice constant by using the lattice dispersion relation. The corre- 
sponding relation for the two-particle levels reads 

(2sinh^^) 2 = r^ + £\ (6.17) 

Since the absolute values of the momenta k v are smaller than the energy levels W v below 
the inelastic threshold (especially for the lowest level with k ~ 0) we can assume that 
the dominating part of the lattice effects is removed by replacing 

— -» 2 sinh ^ ' . (6.18) 
2 2 

As expected the lattice effects grow with increasing distance from the critical point 
(k = k c = 0.30423(1), J — 0). For the lowest energy levels they remain smaller than 
one standard deviation. For the higher energy levels below the inelastic threshold the 
lattice effects are at most twice as large. As we will see in sect. 9 at the simulation point 
(k = 0.320, J — 0.02) with our largest pion mass (and therefore the highest energy values) 
this will give a visible difference in the calculated scattering phase shifts. Compared 
with the lattice effects the influence of the other uncertainties mentioned in this section 
is negligible (see table 6). 

Due to the loss of rotational invariance on the lattice it is not completely unique 
how to calculate k v on the basis of eq. (6.17). We circumvent this problem by using the 
direction of the momentum corresponding to the nearest free energy level jtt(L/2ir) 2 k 2 
listed in table 5. The associated ambiguity might lead to a systematic error which is of 
the same size as the statistical error. 
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6.4 A self-adjusting exponential fit 

We have developed a method to determine the upper end t max of the range of the 
exponential fit (6.16) automatically. The eigenvalues X v (t,t ) are sorted such that the 
corresponding eigenvectors u v (t) of the correlation function matrix V(t,t ) at successive 
values of t are as parallel as possible: Ideally the eigenvector corresponding to an 
energy level W v should be time-independent (see sect. 6.2). To get the optimal assign- 
ment of (almost parallel) eigenvectors is a standard problem of linear optimization. It 
can be solved starting from the matrix A Vfl = (u v (t), u^t+l) ) of the scalar products of 
the eigenvectors on successive time slices: 

max Yl x i™ A w (6.19) 
with Yl x i™ = ^2 X ^ = 1 and x v e i ' 1} • 

v n 

We have approximated the exact solution by the following procedure: We determine the 
largest matrix element and assign u v (t) to u^(i+l). After deleting the corresponding 
row and column we continue in this way with the truncated matrix. 

The sorting of the eigenvalues according to the collinearity of their eigenvectors should 
be equivalent to a sorting with respect to their size, provided the eigenvalues at to+1 are 
arranged in accordance with this criterion. But due to statistical fluctuations at some 
time slice t = t max +l this fails, thus determining the upper end t max of the fit range (for 
an example see fig. 7). 

6.5 Final choice of the method 

Table 6 shows the weak dependence of the spectrum on the choice of variations displayed 
in fig. 6. The variations are done independently for each degree of freedom starting at 
a standard choice. From table 6 we learn to keep the reference time to small and to 
guarantee a complete set of linear independent eigenvectors corresponding to the energy 
levels in the elastic region (f >z). 

Among all the possibilities we selected the following list of parameters, which are 
labeled by a • in table 6: 

• reference time: to — 

• block length: 2048 (with 18 — 40 blocks we obtain enough statistics) 

• choice of the strategy (3) (taking the mean of the energy values) 

• dimension of V(t,t ): at least f = i+ 1 

• temporal lattice extent T = 32 and 40, respectively (24 3 • 60 lattice only with k = 
0.315 and J = 0.01) 

• "Liischer's" two-particle wave function (6.4) with isospin and plane waves (6.3) 
with isospin 1 and 2 

For comparison we also illustrate the systematic shift of the energy values due to the 
replacement (6.18). The following section summarizes the results of our analysis of the 
two-particle energy spectrum obtained along these lines. 
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7 Numerical results for the two-particle energy spec- 
trum 

The correlation function matrices are analysed according to the criteria mentioned above. 
The results are shown in tables 7 to 10 without corrections for lattice effects. These final 
results are the starting point for the calculation of the scattering phase shifts and the 
determination of the resonance parameters m a and T a (and the coupling g^, respectively) 
as discussed in the following section. For comparison we also show the perturbatively 
calculated energy spectrum using as input parameters the results of sect. 9 (see table 11) 
instead of the estimates in table 4. 

Fig. 10 shows all data of the two-particle energy spectra for isospin and 2 at a 
glance, whereas fig. 8 and 9 give an enlarged overview at the simulation point (k = 0.315, 
J = 0.010): 

• Crosses symbolize the measured energy values below the inelastic threshold 4??^. 
The errors are smaller than the symbols and lattice effects are taken into account 
by means of the replacement (6.18). 

• The solid lines passing through our data represent the energy spectrum calculated 
perturbatively on the basis of the results of sect. 9. We observe good agreement 
with perturbation theory. 

• The perturbative spectrum computed from the estimates of table 4 is indicated by 
the dashed curves in fig. 8 

For isospin 1 the picture would look similar to the isospin 2 case except that there is no 
energy level around 2m,r. 

In the isospin-0 resonance channel we see the expected trend (cf. fig. 2): The 
first and second level are coming close to each other displaying the so-called "avoided 
level crossing" such that we get an impression where the resonance plateau might be 
(the dotted line at about W « 3m T ). 

In fig. 8 we also compare our results with perturbative calculations based on the 
estimates of table 4. The reason why the corresponding dashed lines do not fit to the 
data is related to the relatively large width of the resonance (see sect. 9): the naive fit 
to the propagator in momentum space does not yield the correct value for the resonance 
mass. However, there is very nice agreement between perturbation theory and our data, 
if we take the resonance parameters m a and Y a (g#) determined in sect. 9. 

In the isospin- 1 and the isospin-2 channel the numerically calculated two-particle 
energy spectrum agrees very well with the perturbative prediction, which is close to the 
free energy spectrum and depends only weakly on m a (see e.g. fig. 9 for isospin 2). As 
shown in fig. 10 we obtain a similar behaviour at all simulation points. 
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two-particle energy spectrum 

ISOSPIN 





Wo 


W 1 


w 2 


w 3 


W 4 W 5 W 6 




*o = 


0.456(2) 


0.644(2) 0.747(3) 0.902(4) 1.22(2) 


• 


to = 1 


0.452(4) 


0.64(1) 


0.75(1) 


0.86(2) 


1.25(7) 




to = 2 


0.43(3) 


0.60(3) 


0.70(2) 


0.93(5) 


1.8(2) 




hik = 1024 


0.452(3) 


0.643(3) 


0.744(3) 


0.906(4) 


1.22(2) 




hik = 2048 


0.456(2) 


0.644(2) 


0.747(3) 


0.902(4) 


1.22(2) 


• 


hik = 4096 


0.458(3) 


0.646(2) 


0.749(3) 


0.900(4) 


1.21(2) 




hik = 8192 


0.457(2) 


0.649(3) 


0.749(3) 


0.898(5) 


1.21(2) 




r = 4 


0.455(2) 


0.646(2) 


0.750(3) 


0.991(6) 






f = 5 


0.456(2) 


0.644(2) 


0.747(3) 


0.902(4) 


1.22(2) 


■ 


f = 6 


0.455(2) 


0.644(3) 


0.745(3) 


0.896(4) 


1.022(4) 1.44(4) 




f = 7 


0.455(2) 


0.643(3) 


0.744(4) 


0.895(4) 


1.019(4) 1.141(5) 1.50(5) 




® 


0.459(1) 


0.650(3) 


0.748(4) 


0.899(4) 


1.200(9) 






0.455(1) 


0.645(2) 


0.744(2) 


0.899(3) 


1.22(2) 




® 


0.456(2) 


0.644(2) 


0.747(3) 


0.902(4) 


1.22(2) 


• 


24 3 • 32 


0.456(2) 


0.644(2) 


0.747(3) 


0.902(4) 


1.22(2) 


• 


24 3 • 60 


0.454(2) 


0.652(2) 


0.742(3) 


0.909(3) 


1.23(2) 




pi. waves 


0.454(2) 


0.647(2) 


0.751(4) 


0.905(4) 


1.028(4) 




Liischer fct. 


0.456(2) 


0.644(2) 


0.747(3) 


0.902(4) 


1.22(2) 


• 


W v -> 4sinh^ 


0.457(2) 


0.647(3) 


0.751(3) 


0.909(4) 


1.24(2) 





Table 6: Comparison of the two-particle energy data in the isospin-0 channel at (k = 0.315, 
J — 0.01) for different methods of analysis: different reference times t , block lengths hik and 
dimensions f of the correlation function matrix. Further variations concern the choice of the 
wave function and the strategy for calculating the block mean values (fig. 6). Additionally we 
compare lattices with different time extent: 24 3 • 32 and 24 3 • 60 lattices. The lines with • refer 
to our choice of final results (they are displayed several times). 
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two-particle energy spectrum 
isospin 







w 


W 1 


w 2 


w 3 


W 4 W 5 W 6 


16 3 


32 


0.352(2) 


0.548(2) 


1.032(8) 










0.363 


0.548 


0.903 






20 3 


32 


0.358(2) 


0.526(3) 


0.753(2) 


1.15(2) 








0.368 


0.531 


0.759 


0.99 




24 3 


32 


0.360(3) 


0.520(3) 


0.661(3) 


0.849(4) 1.28(5) 






0.371 


0.519 


0.669 


0.853 


0.99 


32 3 


40 


0.368(2) 


0.493(2) 


0.566(3) 


0.681(2) 


0.780(5) 1.12(5) 






0.373 


0.498 


0.570 


0.689 


0.785 0.87 



isospin 2 







W 


W 1 


w 2 


w 3 


W 4 


w 5 


16 3 


32 


0.368(1) 


0.871(1) 1.165(2) 












0.376 


0.877 


1.182 








20 3 


32 


0.369(2) 


0.730(1) 


0.962(1) 


1.142(2) 










0.375 


0.736 


0.971 


1.154 






24 3 


32 


0.363(2) 


0.641(1) 


0.829(1) 


0.973(1) 


1.091(1) 








0.375 


0.646 


0.834 


0.984 


1.114 




32 3 


40 


0.369(1) 


0.538(1) 


0.667(1) 


0.771(1) 


0.861(1) 


0.949(1) 






0.374 


0.544 


0.672 


0.778 


0.871 


0.958 



Table 7: The two-particle energy spectrum for n = 0.310 and J = 0.005 calculated as described 
in sect. 6.5. For comparison we also show the perturbative predictions (numbers without errors) 
based on the results of table 11. Only some higher energy levels above the four particle threshold 
W > 4m x = 0.748(3) (indicated by the zigzag line) show larger deviations. The errors given 
are purely statistical. 



20 




two-particle energy spectrum 
isospin 







W 


W 1 


w 2 


w 3 


W 4 


w 5 


W 6 


16 3 


32 


0.447(3) 


0.703(4) 


0.948(4) 1.36(2) 










0.454 


0.704 


0.953 


1.24 






20 3 


32 


0.454(3) 


0.673(3) 


0.814(4) 


1.024(4) 


1.41(4) 








0.457 


0.678 


0.822 


1.035 


1.20 




24 3 


32 


0.456(2) 


0.644(2) 


0.747(3) 


0.902(4) 


1.22(2) 








0.458 


0.651 


0.748 


0.904 


1.030 




32 3 


40 


0.450(2) 


0.588(2) 


0.667(2) 


0.753(2) 0.834(2) 0.912(3) 


1.21(2) 






0.459 


0.593 


0.676 


0.760 


0.835 


0.916 


1.01 










isospin 


2 








W 


W 1 


w 2 


w 3 


W 4 


w 5 


16 3 


32 


0.457(3) 


0.905(2) 


1.193(3) 


1.413(3) 










0.461 


0.916 


1.211 


1.441 






20 3 


32 


0.458(2) 


0.776(2) 


0.994(1) 


1.168(3) 


1.306(3) 








0.460 


0.782 


1.006 


1.184 


1.340 




24 3 


32 


0.458(2) 


0.694(1) 


0.865(1) 


1.012(2) 


1.127(2) 


1.240(2) 






0.460 


0.699 


0.875 


1.019 


1.145 


1.263 


32 3 


40 


0.453(1) 


0.602(1) 


0.720(1) 0.817(1) 0.901(1) 


0.987(1) 






0.460 


0.606 


0.723 


0.822 


0.911 


0.994 



Table 8: The two-particle energy spectrum for n = 0.315 and J = 0.01 calculated as described 
in sect. 6.5. For comparison we also show the perturbative predictions (numbers without errors) 
based on the results of table 11. Only some higher energy levels above the four particle threshold 
W > 4m,r = 0.920(3) (indicated by the zigzag line) show larger deviations. The errors given 
are purely statistical. 
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k = 0.315 
J - 0.010 



two-particle energy spectrum 
isospin 1 







W 


W 1 


w 2 W 3 


W 4 


16 3 


32 


0.897(2) 


1.181(3) 


1.401(2) 








0.909 


1.201 


1.435 




20 3 


32 


0.769(1) 


0.988(2) 


1.166(1) 1.303(2) 








0.778 


1.000 


1.181 1.340 




24 3 


32 


0.690(1) 0.863(1) 


1.003(2) 1.121(2) 


1.233(2) 






0.697 


0.871 


1.016 1.144 


1.257 


32 3 


40 


0.599(1) 0.714(1) 0.814(1) 0.901(1) 


0.980(1) 






0.605 


0.721 


0.821 0.910 


0.991 



Table 9: The two-particle energy spectrum for k = 0.315, J = 0.01 and isospin 1 calculated 
as described in sect. 6.5. For comparison we also show the perturbative predictions (numbers 
without errors) based on the results of table 11. Only some higher energy levels above the four 
particle threshold W > Am v = 0.920(3) (indicated by the zigzag line) show larger deviations. 
The errors given are purely statistical. 
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k = 0.320 
J - 0.020 



two-particle energy spectrum 
isospin 







w 


W x 


w 2 


w 3 


W 4 


w 5 


W 6 


12 3 


32 


0.575(2) 


0.888(4) 


1.240(3) 1.80(2) 












0.579 


0.894 


1.261 


1.652 








16 3 


32 


0.582(4) 


0.844(5) 


1.031(5) 


1.29(1) 


1.79(6) 










0.585 


0.852 


1.038 


1.30 


1.50 






20 3 


32 


0.582(3) 


0.802(4) 


0.920(4) 


1.096(5) 


1.248(6) 


1.60(6) 








0.587 


0.810 


0.928 


1.109 


1.256 


1.40 




24 3 


40 


0.582(3) 


0.756(3) 


0.861(3) 


0.977(4) 1.101(8) 


1.192(5) 


1.51(6) 






0.588 


0.767 


0.869 


0.992 


1.100 


1.209 


1.33 



isospin 2 







W 


W x 


w 2 


w 3 


W 4 


w 5 


12 3 


32 


0.587(2) 


1.189(1) 1.558(2) 


1.842(3) 










0.591 


1.211 


1.607 


1.914 






16 3 


32 


0.589(2) 


0.978(1) 


1.246(2) 


1.449(2) 


1.608(3) 








0.590 


0.986 


1.265 


1.487 


1.680 




20 3 


32 


0.587(3) 


0.857(2) 


1.059(2) 


1.221(2) 


1.353(2) 


1.487(3) 






0.590 


0.864 


1.071 


1.240 


1.390 


1.529 


24 3 


40 


0.586(2) 


0.786(2) 


0.940(1) 


1.074(2) 


1.185(2) 


1.292(2) 






0.590 


0.790 


0.949 


1.083 


1.203 


1.315 



Table 10: The two-particle energy spectrum for k = 0.320 and J = 0.020 calculated as 
described in sect. 6.5. For comparison we also show the perturbative predictions (numbers 
without errors) based on the results of table 11. Only some higher energy levels above the four 
particle threshold W > 4m x = 1.179(3) (indicated by the zigzag line) show larger deviations. 
The errors given are purely statistical. 
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8 The properties of the wave functions 

For both kinds of wave functions (6.3) and (6.4) we obtain the same energy spectrum 
within a standard deviation. This fact shows how reliable our diagonalization procedure 
is. However, the question remains, if Liischer's wave functions (6.4) have an improved 
projection on energy eigenstates compared with the plane waves (6.3). 

In order to get an idea of the projection properties we consider the normalized spectral 
amplitudes (cf. eq. (6.15)). In figs. 11 and 12 we plot \c v -\ 2 for the isospin-0 operators 
Oj,j = a, 0, 1, . . . , f — 2, with plane waves and Liischer's wave functions, respectively, 
at (k = 0.315, J = 0.01, L = 24). Each box shows \d<\ 2 for a fixed value of j. Within 
statistical errors, the results for \c v a \ 2 are independent of the type of wave function used 
in the operators O ,O 1 , . . . Of course, this has to be so, since O a is the same operator 
in both cases. This agreement provides another consistency check of our calculations. 
In the case of the plane waves we see clearly why it was vital to include the cr-operator: 
O ,O 2 , 3 and 4 generate essentially only the energy eigenstates \v) with v = 0, 3, 4, 5. 
So Oi needs a partner, if |1) is to be separated from |2). 

The operators Oj constructed with Liischer's wave functions do what they are sup- 
posed to do: Oj generates predominantly the energy eigenstate \ j), except that 0\ leads 
to a strong v = 2 component also. This could be due to an inaccurate estimate of qf, 
which had to be fixed beforehand (cf. sect. 6.1). In summary we have to conclude that 
there is no essential difference between Liischer's wave functions and plane waves. 
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9 Numerical results for the scattering phases 

For each lattice extent L and each two-pion energy level W v < 4m v we get one value of 
the scattering phase shift 6$, 6\, and 6%, respectively, in the elastic region < k^/m^ < y/3 
using the key relation (2.4) 

8l{k v ) = -d> ( ^ ) modTT. (9.1) 

The momentum k v corresponding to W v has to be calculated with the help of the energy 
momentum relation (6.17) as discussed in sect. 6.3. 

Fig. 13 illustrates the calculation of the scattering phase shifts for the simulation point 
(k = 0.315, J = 0.01). For some fixed values oim v L we plot the function — (f> (§f ) mod7r, 
already shown in fig. 1, versus k/m^. The value of this function at a momentum k v (W v ) < 
x/3 m v gives the scattering phase shift at that momentum. The data points in fig. 13 
show the momenta extracted from the two-pion energy levels and the corresponding 
phase shifts whose errors are determined by the slope of <j). Each passing of the function 

— (f) through a multiple of 7r turns out to correspond to one energy level W v . 

Comparing the different values of L we observe the following: As L increases, more 
and more energy levels lie below the inelastic threshold. At the same time the function 

— (f) (§7) becomes steeper when considered as function of k. This fact has important 
consequences for the error propagation: On larger lattices higher accuracy of the energy 
values W v is required, though it is increasingly harder to achieve. 

There is another consequence of the steepness of the function —(f). The difference 
between the energy values determined with or without lattice effects taken into account 
is in general not much larger than the statistical error. Especially at the simulation point 
(k = 0.320, J = 0.02), where the largest energy values appear, this difference leads to a 
systematic shift of the resulting scattering phases due to the large slope of the function 

— (f). Fig. 14 illustrates this effect: Using the lattice dispersion relation (6.17) instead of 
the continuum dispersion relation the scattering phases are shifted downwards along the 
dotted lines representing —(f) mod7r. 

How do our results compare with perturbation theory? In fig. 15 we plot all values 
for 6® that we obtained at (k = 0.315, J = 0.01) together with the perturbative prediction 
based on the estimates m a and c/r (dashed curve). The comparison shows that the value 
of the resonance mass m a lies below the estimate m a from the fit to the er-propagator in 
momentum space. 

Figs. 17 and 16 summarize our results for all isospin values and simulation points. For 
isospin 1 and 2 the agreement of the numerical results and the perturbative prediction 
is quite good (even if we use the estimates m a and of table 4). In these cases the use 
of the lattice dispersion relation (6.17) is essential. Note that the data fall onto a single 
curve, although they originate from different lattice sizes. This justifies a posteriori the 
assumption that the influence of higher angular momenta and finite volume polarization 
effects can be neglected. 



25 



Finally we discuss different methods for the determination of the resonance mass m„ 
and the decay width r CT from the measured scattering phases: 

1. Breit-Wigner-fit (model independent): Fitting the values of isospin-0 scattering 

phases around the resonance (7r/4 < 8$ < Stv/4) by means of the (relativistic) Breit- 
Wigner formula 

tm ( 4 °-2) = ^rr <9 - 2) 

it is possible to calculate the resonance mass m a and its width r CT . One major 
advantage of this fit is the determination of the resonance parameters m a and r CT 
without additional assumptions concerning the model under investigation. The fit 
works reasonably well although there are only few data points inside the resonance 
region (vr/4 < <S ° < 3tt/4). 

2. Perturbative fits (model dependent): Outside the resonance region the Breit- 

Wigner-formula (9.2) does not apply. Using the perturbative formula (3.6) 

S°o - Si + 81 (9.3) 



J 0,r 

... , co N ~ 1 ml -ml k 

with tan« 0iJ =^— — ^- 



Qr ml - ml , / 4/c 2 + ml \ N + 1 k 



and C = ^) 



96vr fcVF V / aa 48tt VF 

it is possible to fit all data points in the elastic region by a two-parameter fit with 
respect to the resonance mass m a and the coupling g^. 

The width r CT of the er-resonance may then be calculated from the perturbative 
formula [211: 



9r- 



N - 1 ml -ml 



yjml - Ami ■ ( 9 - 4 ) 



Qr = ZZ y *> . (9.5) 



967T ml 

If in addition to the pion mass m v also the values of the wave function renorma- 
lization constant Z and the field expectation value S are taken to be known from 
table 4, the same ansatz (9.3) can be used in order to perform a 1-parameter fit 
with respect to m a provided that we replace by its definition 

Xml_-ml) 

The results of the latter fit serve as basis for the perturbative calculations of the 
energy spectrum (see tables 7 to 10) and the computation of the curves shown in 
our figures. 

The results for the resonance mass m a , the decay width r CT and the coupling obtained 
by means of these fits are summarized in table 11. If one of the variables gR or r CT is not 
a fit parameter it is calculated by means of the formula? (9.5) and (9.4), respectively. 
The values do not depend very much on the kind of fit ansatz used. The results for the 
decay width agree reasonably well with the perturbative predictions, while the resonance 
masses m a lie systematically (about 5 %) below the estimates m a from the fit to the 
propagator in momentum space. However, this discrepancy should not be too surprising, 
because the form of the propagator used in the determination of m a corresponds to a 
stable particle and the width of the er-resonance for our choice of parameters is rather 
large: r ff sa 0.11-0.18 m a . 
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SIMULATION POINTS 



K 
J 


1 

0.310 
0.005 


2 

0.315 
0.010 


3 

0.320 
0.020 


Estimates 


rh a 

t a (m a ,m v ,g R ) 


0.538(1) 
22.0(1) 
0.074(1) 


0.720(1) 
24.4(1) 
0.121(1) 


0.906(2) 
27.0(1) 
0.165(1) 


Breit-Wigner-fit 


m a 

g R (m a ,m^,T,,Z) 


0.520(1) 
20.3(1) 
0.079(4) 


0.706(2) 
23.4(1) 
0.130(9) 


0.885(2) 
25.8(2) 
0.160(17) 


Perturbative fit to 


& 9r 






m a 
9r 

T a (m a ,m v ,g R ) 


0.513(1) 
21.1(8) 
0.064(3) 


0.688(2) 
25.7(1.0) 
0.116(5) 


0.863(3) 
28.1(1.4) 
0.156(8) 


Perturbative fit to 


m a 

g R {m a ,mT !) Y>,Z) 
T a (m a ,m 7r ,g R ) 


0.515(1) 
19.8(1) 
0.060(1) 


0.691(1) 
22.3(1) 
0.102(1) 


0.868(2) 
24.6(1) 
0.138(1) 



Table 11: Results for the resonance mass m CT , the coupling g R and the decay width T CT of 
the (T-particle obtained by fits with the relativistic Breit-Wigner-formula (9.2) and by one- or 
two-parameter fits with the perturbative ansatz (3.6), respectively. If the coupling g R or the 
decay width r CT are not fit parameters, they are calculated by means of eqs. (9.5) and (9.4), 
respectively. These results are compared with the estimates of table 4. The errors given are 
purely statistical. 
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10 Conclusion 

The main result of our paper is that Liischer's method for studying particle scattering 
in massive quantum field theories on finite lattices can be applied successfully even in 
four dimensions. It is numerically very stable, i.e. the results vary only within statistical 
errors under a change of the parameters of the analysis (see sect. 6.3) and it passes all 
our selfconsistency checks (see e.g. sect. 8). Furthermore it was confirmed a posteriori 
that correction terms to the key relation (2.4) resulting from polarisation effects and con- 
tributions of higher angular momenta can be neglected in our simulations (see sect. 7). 
Scaling violations seem to be compensated (at least partially) by using the lattice di- 
spersion relation instead of the continuum dispersion relation, as already observed for 
single-particle states. 

In the isospin-0 channel the expected cr-resonance is found. Its mass and width can 
be extracted with reasonable accuracy. Measuring the mass by a fit to the propagator in 
momentum space gives - not surprisingly - a slightly different value. The discrepancy 
is, however, small enough to leave the previous calculations of an upper bound to the 
Higgs boson mass essentially unchanged (unless the difference increases considerably in 
the limit J — > 0). 

The other channels show no resonance structure. In accordance with [28] we see no 
p-resonance. Hence the formation of the p-meson cannot be understood in the framework 
of the $ 4 model. However, this conclusion has to be considered with some caution. If 
the a--operator is omitted in the study of the isospin-0 channel, all signs of a resonance 
disappear. Could it be that we have left out the "essential" operator in the isospin-1 
channel? 

Since one might argue that an operator which is more strongly localized than the 
plane waves used above should have an enhanced projection on the resonance, we have 
added in the analysis of the isospin-1 channel the operator 



which corresponds to the conserved 0(3)-current in the continuum limit. But its inclu- 
sion does not alter our results. 

As already observed in other investigations, due to the triviality of the $ 4 -theory in 
four dimensions renormalized perturbation theory is reliable. We find complete agree- 
ment with our nonperturbative results, provided the renormalized parameters, especially 
m a , are chosen properly. This finding corroborates our above statement that Liischer's 
method works well. Furthermore, in view of the agreement with perturbation theory, it 
seems unlikely to us that we have overlooked a resonance in the isospin-1 channel, which 
would be a nonperturbative phenomenon. 

The successful application of Liischer's method to a four-dimensional quantum field 
theory as well as the above mentioned investigations of lower dimensional systems give 
rise to some optimism concerning the study of more demanding models, in particular 




L 




QCD. 
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